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Abstract 

We carry out Monte- Carlo simulations to discuss critical properties of a classical two- 
^ , dimensional XY frustrated helimagnet on a square lattice. We find two successive phase 

transitions upon the temperature decreasing: the first one is associated with breaking of 
a discrete Z2 symmetry and the second one is of the Berezinskii-Kosterlitz-Thouless (BKT) 
type at which the SO{2) symmetry breaks. Thus, a narrow region exists on the phase diagram 
S ' between lines of the Ising and the BKT transitions that corresponds to a chiral spin liquid. 
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5 : 1. Introduction 

Frustrated magnets have attracted much attention in recent years. Exotic spin hquid phases are 
of special interest which have been found in some of them jl . A chiral spin liquid phase is an 
example of such an exotic state of matter in vi^hich there are neither quasi long-range nor long-range 
Xjrs , magnetic orders but a chiral order parameter (S^ x Sj) is nonzero. Existence of such a phase is 

C^ ' discussed in context of one-dimensional frustrated quantum magnetic systems [5], and it is found 

Cn . experimentally in Ref. [3]. 

^^ ' In larger dimensions, one of the systems in which the chiral spin liquid phase can be found 

CO . at finite temperature is a classical planar (XY) helimagnet with Z2 5*0(2) symmetry in which 

^^ ' the helical structure results from a competition of exchange interactions between localized spins. 

^ ^ . Critical behavior of spin systems from this class is described by two order parameters. Besides 

the conventional magnetization with 50(2) symmetry, one has to take into account also the chiral 

order parameter that is an Ising variable with Z2 symmetry. This parameter characterizes the 

direction of the helix twist and distinguishes left-handed and right-handed helical structures. 

%^ , In three-dimensional helimagnets, the phase transitions on the magnetic and the chiral order 

^.' parameters occur simultaneously. It was found numerically that the transition is of the weak first 

order or of the " almost second order" Hj [5] type in helical antiferromagnets on a body-centered 

tetragonal lattice [6] and on a simple cubic lattice with an extra competing exchange coupling 

along one axis [?]■ These systems belong to the same (pseudo) universality class as, e.g., the 

model on a stacked-triangular lattice i8j and V2,2 Stiefel model [9.. The possibility of existence 

and stabilization of the chiral spin liquid phase by, e.g., Dzyaloshinsky-Moria interaction in 3D 

helimagnets is discussed recently in Refs. '10) . 

In two dimensions, the situation is rather different |11) . Two successive transitions were ob- 
served with the temperature decreasing. The chiral order appears as a result of the first transition 
that is of the Ising type. Another one is the Berezinskii-Kosterlitz-Thouless (BKT) transition 
driven by the unbinding of vortex-antivortex pairs [12, . Then, the chiral spin liquid phase arises 
between these transitions with the chiral order and without a magnetic one. Various 2D systems 
from the class Z2 SO{2) were investigated numerically (see Ref. [T^ for review): triangular 
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antiferromagnet [131 [TS], J1-J2 model [TB], the Coulomb gas system of half-integer charges [T7] . 
two coupled XY models [T^ , Ising-XY model [THl [131 [2D1 and the generalized fully frustrated XY 
model [Hj. And surely, the most famous of them is the fully frustrated XY model (FFXY) in- 
troduced by Villain [35]. This model is of great interest because it describes a superconducting 
array of Josephson junctions under an external transverse magnetic field [23'. It was found that 
the temperature of the Ising transition T/ is 1-3% larger than that of the BKT transition for the 
most of above-named systems [23l [TTl [231 [13] . 

Korshunov argued [25j that a phase transition, driven by unbinding of kink-antikink pairs on 
the domain walls associated with the Z2 symmetry, can take place in models similar to 2D FFXY 
one at temperatures appreciably smaller than Tbkt (see also Ref. [26]). Such a transition could 
lead to a decoupling of phase coherence across domain boundaries, producing in this way two 
separate bulk transitions with Tbkt < Tj [57]. It was pointed out however in Ref. [25 that these 
two continuous transitions can merge into a single first order one. These conclusions do not depend 
on the particular form of interactions in system as soon as the ground state degeneracy remains 
the same. They are confirmed by numerical studies of the models mentioned above [23l [TTl [131 fT3l 
[16l[ll[15l[T7]. 

Nevertheless the situation remains contradictory in 2D helimagnets belonging to the same 
Z2 18) SO{2) class as FFXY model and the antiferromagnet on the triangular lattice. Garel and 
Doniach [55] (see also [5S]) considered the simplest helimagnet on a square lattice with an extra 
competing exchange coupling along one axis that is describing by the Hamiltonian 



^ = E 



(jl COs((^x - (/Jx+a) + J2 COs((px - <(5x+2a) - Jb COs((^x - V'x+b) j , (1) 



where the sum runs over sites x — {xa,Xb) of the lattice, a — (1,0) and b = (0,1) are unit 
vectors of the lattice, the coupling constants Ji_2 are positive. Using arguments of Ref. [30^, they 
concluded [28] that at low temperatures the vertices are bound by strings, which would inhibit 
the BKT transition and make the Ising transition occur first with the temperature increasing. 
Kolezhuk noticed [31^ that those arguments are not valid for a helimagnet, and showed that the 
Ising transition temperature is larger than the BKT one at least near the Lifshitz point J2 = Ji/4. 
It was found by Monte Carlo simulations in the recent paper [35' that Tbkt > Tj at J2 = 0.3 
and Jl = Jb = 1 (i.e., very near the Lifshitz point) in accordance with Ref. |i28j and in contrast to 
Ref. [S]. 

To account for the discordance between results for helical magnets and the general arguments 
for Z2 (E) SO{2) class, we perform extensive Monte Carlo simulations of the model (1) for different 
values of J2. We obtain reliable results at J2 > 0.4Ji which show that Tbkt < Tj. On the other 
hand the value of Tj close to the Lifshitz point is hiding among effects of the finite size scaling and 
is not accessible for ordinary estimation methods. We obtain the Ising transition temperature from 
the chiral order parameter distribution and find that Tbkt < Tj near the Lifshitz point too. At the 
same time we find in accordance with results of Ref. [32' that the specific heat and susceptibilities 
have subsidiary peaks at low T < Tbkt near the Lifshitz point. These are anomalies which are 
attributed in Ref. [32] to the Ising phase transition. However, we demonstrate that these anomalies 
do not signify a continuous phase transition. Apparently, their origin is in metastable states which 
lead also to a peculiar distribution of the chiral order parameter. We find no such features in the 
specific heat and susceptibilities far from the Lifshitz point (at J2 > 0.4Ji). As a result we obtain 
the phase diagram shown in Fig. [T] 

The rest of the present paper is organized as follows. We discuss in Sec. [5T] the model ([1]) 
in more detail and introduce quantities to be found in our calculations. Numerical results are 
discussed in Sec. [33 In particular, the Ising and the BKT transitions are considered in Sees. 13.1.1 
and [323 respectively. The neighborhood of the Lifshitz point and the phase diagram are discussed 
in Sec. 13.3.1 Sec. |3D contains our conclusions. 
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Figure 1: Phase diagram of the model ^ that is found in the present paper. 



2. Model and methods 

We consider the model (H]) of the classical XY magnet on a square lattice. We set Ji ~ Jb = 1 
for simplicity and the value of the extra exchange interaction J2 is a variable. The Lifshitz point 
corresponds to J2 = 1/4 in this notation. The system has a collinear antiferromagnetic ground 
state at J2 < 1/4. To discuss the phase transition from the (quasi-) antiferromagnetic phase to the 
paramagnetic one we consider J2 = and J2 = 0.1 (see Fig. [ij. The ground state has a helical 
ordering at J2 > 1/4. The turn angle Oq between two neighboring spins along a axis is given by 
cos 00 = —Ji/4:J2 at zero temperature. 

To discuss the number and the sequence of phase transitions from the (quasi-)helical phase to 
the paramagnetic one we consider J2 « 0.309, 0.5 and 1.76 corresponding at T = to angles of 
commensurate helices Oa = 47r/5, 2tt/3 and Gvr/ll, respectively. We use lattices with L^ cites, where 
L is divisible by the size of the helix pitch and it lies in the range from 20 to 120. We apply the 
periodic (toric) boundary conditions as well as the cylindrical ones (i.e., with the periodic condition 
along the b axis and the free one along the a axis). We have found that both conditions lead to 
the same values of transition temperatures and indexes. In contrast values of Binder's cumulants 
and the chiral order parameter distribution at J2 ~ 0.309 depend on boundary conditions as we 
discuss below in detail. Standard Metropolis algorithm [35] has been used. The thermalization 
was maintained within 4 • 10^ Monte Carlo steps in each simulation. Averages have been calculated 
within 3.6 • 10^ steps for ordinary points and 6 • 10® for points close to the critical ones. We have 
used also the histogram analysis technique in which the range of each quantity has been divided 
into 6.4 • 10^ bins. 



2.1. Order parameters 

The BKT transition is driven by the magnetic order parameter for which we use two definitions. 
Similar to the triangular lattice [14] , one can introduce a number of sublattices in the case of helix 
pitches which are divisible by the lattice constant. Then, one can write for the magnetic order 
parameter 






■Ek^)' 



(2) 



where index i enumerates Usi sublattices, the sum over x^ runs over sites of the i-th sublattice, 
spin Sxi = (cos (/fxi J sin 0xi ) is a classical two-component unit vector, and (. . .) denotes the ther- 
mal average. The second definition of the order parameter is valid both for commensurate and 



incommensurate helices 

M, = Tl.S,a+,b, M= /-V(M^). (3) 




Our calculations show that definitions ([2]) and ([3|) lead to the same results away from the Lifshitz 
point. We have found that m shows an anomalous behavior at J2 « 1/4 and we use definition 
(J3|) in this case. Thus, we demonstrate that it is useful in numerical discussion of helimagnets to 
choose parameters of the Hamiltonian so that the helix pitch at T = to be commensurate. 
The Ising transition is driven by the chiral order parameter defined as 

^= r2 ■ a y!sin((^x-'^x+a), k ^ ^JW) . (4) 

2.2. Susceptibilities and cumulants 

We introduce corresponding susceptibilities for all order parameters [36] 

L2 



i.((P^)-(bl>^)> 



T < n. 



Xp={ ^ (5) 

The second line in this definition is used below for estimation of critical exponents. Binder's 
cumulants |36| are define as 

3(p2 
We discuss also the cumulant 



"r=i- iB-2^ C) 



using which the critical exponent v^ can be found by finite-size scaling analysis [3 7) . 

2.3. Helicity modulus 

It is useful to introduce the helicity modulus (or the spin stiffness) [38] to discuss the BKT transition 
that is defined by the increase in the free energy density F due to a small twist A^ across the 
system in one direction (a or b) 



' dA, 



(8) 



where ^ = a,b denotes the direction. Important universal properties of a BKT transition predicted 
by Kosterlitz and Nelson 39, are the jump of the helicity modulus ^ from zero at T > Tqkt to 
the value of ^Tbkt/t^ at T = Tbkt and the value of the exponent ri{T = Tbkt) = 1/4. These 
properties have become standard methods of finding the transition temperature. 

As a result of the fact that the exchange couplings along a and b axes are different, the helicity 
moduli in these directions differ too. Thus, at zero temperature Ta(0) — AJ2 — Ji/(4J2), while 
Tb(0) = Jb- Nevertheless, both Tq and Tb must vanish at the same temperature with the identical 
value of the jump. 

One finds after trivial calculations using Eqs. ^ and (|8]) that the helicity modulus Tf, is 
expressed via correlation functions and has a common view [40) 
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Figure 2: Distribution of the value E',, defined in Figure 3: Distribution of the value E'^ defined in 
Eq. @ for Ja = 0.5, L = and three T values: T > Eq. ^ for Ja = 0.5, T = 0.67 < Ti and different 
Ti, T < Tbkt, and Tbkt <T <Ti. L. 



where T is measured in units of Ji = Jb = 1, we set ks = i, E'^^ = L ^X]x^^'^('^x — Vx+h) and 
X^x cos((px ~ tpx+b)- Similar calculations give for the helicity modulus in the a direction 



K' = L-' 



Ta^{E':)-^{iE',r) + ^{E',) 



(10) 



where E'^ = L ^ Y^^i^^^if^ - V'x+a) + 2 J2 sm{ip^ - (y3x+2a)) and E'^ = -L ^ J2^{cos{(p^ - tpx+a) + 
4 J2 cos(v3x ~ </5x+2a))- It may seem that the last term in Eq. (TTOl) can be discarded as it is done with 



the corresponding term in Eq. (JH) which is equal to zero. However it is not so because (£'^) = at 
all T whereas {E'^) = at T > T/ and {E'^) =^ at T < T/. To demonstrate this let us apply an 
infinitesimal twist Aa across the system in the a direction, i.e., let us replace in Eq. ^ (p^ — </?x+a 
by fx — Vx-i-a + Aa- One writes in the first order in A^ 



y^ (cOs((y9x - (/3x+a + Aq) + J2 COs((^x " '/3x+2a + 2Aa) 

X 

~ ^(cOs((px - <i2x+a) + J2 COs((px " </'x+2a) 

X 

- Aa ^ (sin((^x - '/'x+a) + 2 J2 sin((^x - '/'x+ 



(11) 



Comparing the last term in Eq. (|lip with the chiral order parameter k definition ^ and noting 
that one can use an equivalent definition k — L^^ ^^ sin((^x ~ VxL+2a) we conclude that the last 
term in Eq. (|lip is a linear combination of k and k. However, k and k have opposite signs in the 
case considered. In particular, their combination k + 2J2k in the last term in Eq. (jlip is equal to 
zero at T = 0. Nevertheless our numerical results presented below show that this combination is 
not equal to zero at T 7^ and it can be considered as the Ising order parameter aXT ^ Tj. Then, 
one see from Eq. (fTTj) that A^ plays the role of the "chiral" field and, consequently, dF/d/S.a\Aa=o 
(that is equal in our notation to {E'^)) is proportional to the chiral order parameter which is equal 
to zero at T > Tf and which is finite at T < T/. 

To illustrate this consideration we present in Fig. [2] the distribution of i?^ for J2 = 0.5 and 
for various temperatures below Tbkt ~ 0.671, between Tbkt and Tj w 0.69, and above Tj (the 
values of T/ and Tbkt are obtained below). The distribution has a Gaussian form with the zero 
expected value {E'^} = 0. Fig. [3] shows the distribution of E'^ for T = 0.67 < T/ and various lattice 
sizes. A non-zero expected value of {E'^) is seen. One can observe a double-peak structure of E'^ 
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Figure 4: Same as in Fig. |3] for T >Ti. 



distribution (see Fig. [3]) for small lattices or at T that is close enough to T/, when the system can 
tunnel to a configuration with opposite chirality. The probability of such tunneling is estimated as 



p(/^^ 



k-) 



exp 






(12) 



where fdw is a domain wall tension |41j that is positive at T < Tj and it vanishes at T = T/. That 
is why we observe two peaks for small lattices and only one peak for large ones. Quite expectedly, 
we find the single-peak distribution of E'^ a.t T > Tj demonstrated in Fig. SI 

It should be stressed that the disappearance of the double-peak structure at the critical tem- 
perature Ti{L) is a signature of the transition on a lattice with size L. The value of Tj{L) is 
close to the correct value of the transition temperature for large lattices. We use this circumstance 
below in our analysis of the Lifshitz point neighborhood. 

It should be noted also that we replace in our numerical calculations {E'^) by (Ifi'al) at T < Tj 
in the last term in Eq. (jlOp as it is usually done in considerations of order parameters [36j . It is 
done because the order parameter distribution has tails in both positive and negative regions even 
below the transition temperature. The value (p) is replaced by {\p\) in Eq. ([5]) by the same reason. 

3. Numerical results 

We discuss in this section in detail our results for the special case of J2 = 0.5 that corresponds at 
T = to 120° helical structure with three sublattices. Then, we discuss the phase diagram. We 
consider lattices with L = 24, 30, 36, 42, 48, 60, 72, 90, 120. The lattice with L = 18 is also used 
for estimation of some quantities. 



3.1. Ising transition 

To obtain the Ising transition temperature T/ we use the Binder cumulant crossing method [36] . 
We find for L = 24 and 30 the temperature T^' as a function of ln~ {V /L) at which curves Uk{L) 
intersect for different lattice sizes L' > L. Extrapolation to the thermodynamic limit L' — > 00 
gives the following transition temperature (see Figs. [5] and [6]) 

0.689(1). 



Tr 



(13) 



The dispersion in the value of Tj obtained for the different lattice sizes L gives the error of the 
transition temperature estimation. Notice that error bars arc not shown in Figs. [5] and [6] and in 
all figures below if they are smaller than or comparable with symbols size. 
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Figure 5: Estimation of the transitions tempera- 
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Figure 6: Binder's cumulant Uk{L) defined by 
Eq. ([6| as a function of temperature, for L = 
24, ...,90. 



A peak in the specific heat shown in Fig. [7] is found approximately at the same temperature. 
For the largest lattices the peak is located in the range of temperature from 0.690 to 0.695. One 
expects that this peak corresponds to the logarithmic divergence of the specific heat that is a 
characteristic of the 2D Ising model in which the critical exponent a is equal to zero. Insufficient 
accuracy of our data for the specific heat prevents us from the immediate estimation of a. 

Critical exponents Vk , Pk and 7^ are obtained by the finite-size scaling theory. To estimate the 
exponent Vk , we find a maximum of the quantity Vk given by Eq. ([7|) as a function of lattice size 

L m 



K). 



L^/'^K 



The fitting presented in Fig. [8] gives 



Uk = 0.97(4) 



(14) 



(15) 



that coincides within the computational error with the exact value of :/ = 1 for the 2D Ising model. 
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Figure 7: Specific heat C{T). 



Figure 8: Estimation of the exponent Uk using the 
cumulant 14 defined by Eq. ((7|). 
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Figure 9: Estimation of the exponents l3k and -fk- Figure 10: Binder's cumulant Um{L) as a function 

of temperature, for L = 24, . . . , 90. 



Exponents /3k and 7^ are found from scaling properties of the order parameter k and the 
susceptibiUty Xk at the critical point 



with the following result (see Fig. IH]): 

/3fc= 0.118(8), 7fe = 1.70(6). 



/:^k 



(16) 



(17) 



These values coincide within computational errors with exact values of 1/8 and 7/4, correspond- 
ingly, of the 2D Ising model. Using the scaling relations we find other exponents 



a 



2 - 2uk = 0.06(8), ??fc - 2 - 7fe/i/fe - 0.25(5). 



(18) 



Note that the scaling relation a + 2/3^ + jk ^ 2.00(8) w 2 is satisfied within the computational 
error. 

We have found also the universal value of the Binder cumulant U* — 0.615(6) at the critical 
temperature (see Fig. [S]) that is in agreement with the value U* « 0.611 observed in the 2D Ising 
model [42" with periodic boundary conditions. 

3.2. BKT transition 

According to the Merniin- Wagner theorem |43| there is no spontaneous magnetization at finite 
temperature in 2D magnets with short range interactions and a continuous symmetry. But a 
quasi-long-range order appears at non-zero temperature due to the Berezinskii-Kosterlitz-Thouless 
mechanism [12 in XY magnets with 5*0(2) symmetry. 

It is important that as long as we measure the temperature in units of Ji , the universal value of 
the jump 2Tbkt/t^ perturbs by the factor of Ji/ Jeg, where Jcs is an effective coupling constant. 
The competing exchange coupling J2 gives rise to this factor. To obtain it we considered the 
Coloumb-gas representation of the model ^ using standard duality transformations [12l [44] . As 

a result we obtain 

Jeff = VJb {Ji - 4J2) 



for the antiferromagnetic phase (J2 < Ji/4), and 



J,ff - J^Jb UJ2 - ij^ 



(19) 
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Figure 11: Helicity modulus T;, in the b direction 
and its intersection with the line 2T/-k,JcS. 




Figure 12: Root-mean-square fit error Ac of the 
helicity modulus T(, to the Weber-Minnhagen scal- 
ing equation (pTj) . 
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Figure 13: Same as in Fig. [TT] for Ta 
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for the helimagnetic one (J2 > Ji/4). Eq. (|19p is in accordance with results of Refs. [IHlllS]- In 
particular, Jeff = a/3/2 for J2 = 0.5. 

A few authors have investigated the properties of the Binder cumulant for magnetization as 
an alternative method of the BKT transition temperature estimation [IS] . Due to finite-size cor- 
rections this method gives a value for the transition temperature Tb that is slightly larger than 
the true value Tbkt- Nevertheless this method is useful as it provides an estimation of the BKT 
transition temperature and an extra evidence of separated transitions (if one finds that Tb < Tj). 
Using the Binder cumulant crossing method described above we find Tb = 0.679(2) (see Figs. [S] 
andlini) that is 1.8% smaller than T/ given by Eq. P^ . 

To obtain Tbkt precisely we use the Weber-Minnhagen finite-size-scaling analysis [47 that is 
based on consideration of logarithmic corrections to the value of the helicity modulus at tempera- 
ture close to Tbkt having the form 



T(r,i) 



2T 

TrJcff 



1 



21nL 



(21) 



where c is a fitting parameter. Fixing T we find the root-mean-square error Ac of the least-square 
fit of our numerical data for T(T, L) with different L < 60 that is based on Eq. ((2T|) . The minimum 
of Ac as a function of T gives the value of the transition temperature [47] . We obtain for the helicity 
modulus in the b direction (see Figs. fTTJ and [T2|) 
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Corresponding results for T^ are shown in Figs. [T3] and [TH For the largest lattice of i > 
the Weber-Minnhagen finite-size scaling analysis estimates the BKT transition temperature as 



T; 



BKT 



0.673(2) 



(23) 



that is in agreement with Eq. (|22|) . It should be noted that inaccuracy in estimation of Ta is larger 
than of Tfe. That is why we use below the more precise value (|22|) for comparison between different 
methods. 

To verify our results, we use also cylindric boundary conditions with the periodic condition 
along the b axis and with the free one along the a axis. We obtain results consistent with those for 
periodic conditions. In particular, transitions temperatures Tbkt = 0.671(2) and T/ = 0.6907(6) 
were estimated by the Weber-Minnhagen analysis and the Binder cumulant crossing method (see 
Fig. [T5|) . The universal value of the Binder cumulant at the critical temperature is U* = 0.496(7) 
that is close to the value expected for the 2D Ising model [12] with mixed (cylindric) boundary 
conditions. 

Another indication of the BKT transition is the equality to 0.25 of the exponent rj{T). Below 
the transition temperature the susceptibility diverges with the size of the system as 



xAT,L)^L^-^^^\ 



(24) 



The 



exponent 77 as a function of temperature found using Eq. ()24[) is shown in Fig. 1151 The 

intersection of rj{T) with the bound ry = 0.25 gives 



4kt = 0.676(2). 



(25) 



Comparing Eqs. (|22|) and ((25|) one notes that Tj 



'BKT 



> T, 



(T) 



BKT- 



Because 



v{Tj 



BKT 



Tbkt Jcfi/'^T^Ji^ {Tbkt): we can not exclude that at the true transition temperature the expo- 
nent 77 and the jump of the helicity modulus have non-universal values. Such a possibility has 
been considered for other models from the class Z2 ® S0{2). [331 [Ml IHl HZl SF: If it is so r/ 
has a value smaller than 0.25 and the jump is greater than 2Ji/7rJeff. Thus, our data show that 

v{t1^i^t) « 0.22. 



3.3. Neighborhood of the Lifshitz point and the phase diagram 

Besides the case of J2 = 0.5 considered above in detail we have carried out similar discussions of 
J2 = 0, 0.1, 0.309 and 1.76 to obtain the phase diagram shown in Fig. [H Some results of this 
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'1 '1 

Figure 17: Specific lieat C{T), helicity modulus T;, and susceptibilities Xk.m for J2 « 0.309. 



Table 1: Some results of our discussion of the model ([T]). Here Lmax is the maximum value of L 
considered. 



J2/J1 



Ln 



Tbkt/Ji Tj/Ji 



40 

0.1 30 

: 0.309 47r/5 100 

0.5 27r/3 150 

d 1.76 67r/ll 66 



0.891(2) 
0.781(3) 

0.443(5) 0.48(1) 
0.671(1) 0.690(1) 
1.24(1) 1.285(7) 



consideration are summarized in Table [T] The case of J2 = corresponds to the well known XY 
model on a square lattice and we find Tbkt = 0.891(2) that is consistent with the previous results 

m- 

It should be stressed that we obtain T/ > Tbkt at J2 > 1/4 in accordance with conclusions of 
Refs. [25l[3T] and in contrast to Refs. [281I32]. Because this our finding is at odds with that of the 
similar numerical consideration of the same model carried out in Ref. [32j , our special interest is 
to consider the case of J2 ~ 0.309 that is close to J2 — 0.3 discussed in Ref. [32] . 

The case of J2 ~ 0.309 corresponds at T = to 0o = 47r/5. In particular, we find Tbkt = 
0.443(5) that is very close to the value reported in Ref. [32]. Estimating the temperature of the 
Ising critical point by the Binder cumulant crossing method, we encounter the anomalous behavior 
of the chiral order parameter and do not obtain a reliable result. Apparently, it is the reason why 
the authors of Ref. [32] base their conclusion about the Ising transition on the behavior of the 
specific heat C{T) and susceptibilities Xfe.m- Our results for these quantities are shown in Fig. [TTl 
and they are consistent with those of Ref. [32] • It is seen from Fig. [T7]that the chiral susceptibility 
has a high peak at T « 0.4, while the specific heat and the magnetic susceptibility have subsidiary 
peaks at T « 0.4 which grow with the lattice size increasing. These anomalies at T « 0.4 < Tbkt 
are attributed in Ref. [32] to the Ising transition. 

However, we observe that for other J2 > 1/4 the specific heat has only one peak corresponding 
to the logarithmic divergence which characterizes the Ising transition (see, e.g., Fig.[7]for J2 = 0.5). 
Therefore, the behavior of C{T) shown in Fig. [17] is not normal for the model discussed and it 
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Figure 18: The chiral order parameter distribution 
for J2 = 0.5 and L = 42. 



Figure 19: The chiral order parameter distribution 
for J2 = 0.309, 1/ = 45 at diflerent temperatures. 
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Figure 20: The chiral order parameter distribution tt r,i m, , • n , , ,■ , •, ,■ 

l^igure/1: i he chiral order parameter distribution 
for J, = 0.309, L = 25 and L = 35 at T/J, = 0.43. ^^^ j^ ^ ^^^^^ L = 65 and L = 100 at T/J, = 0.43. 

is characteristic of the Lifshitz point neighborhood. To account for this anomaly we examine the 
behavior of the chiral order parameter in detail. 

Fig. [TH] shows the chiral order parameter distribution for J2 = 0.5 and L = 45. It looks 
like a customary order parameter distribution in a system with the second order transition. In 
particular, the distribution has a Gaussian form below a critical temperature with a peak at k (see 
the curve for T = 0.6). In approaching to the critical temperature the distribution acquires an 
appreciable tail (see curves for T — 0.68 and T = 0.7). Such a broad distribution leads to a peak 
in the susceptibility. The distribution has a peak at fc = above the critical point (see curves for 
T = 0.72 and r = 0.8). 

However the picture for the chiral order parameter distribution is quite different close to the 
Lifshitz point, e.g., at J2 ~ 0.309. We show in Fig. [11] the chiral order parameter distribution for 
L = 45 and different temperatures. One can see one Gaussian peak at T < 0.35 in agreement with 
the common picture described above. But a few additional peaks arise at T > 0.35. Then, the 
number and the breadth of peaks depend on the lattice size L (see Figs. [20] and [21]) and, what 
is much more important, on the boundary conditions. For the cylindrical boundary conditions, 
these peaks are broader and they are accompanied by great number of accessory peaks. Such 
a distribution of the chiral order parameter is characteristic in the case of L == 45 to the range 
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Figure 22: Distribution of E'^ defined in Eq. ([TO 
for L — 65. 



Figure 23: Estimation of Ti for J2 = 0.309 that is 
based on Eq. ((26)) . 



of temperature from T « 0.35 to T « 0.55. One can see from Fig. [19] that at T = 0.57 the 
distribution has a form that is typical for a disordered phase with the peak at zero value of the 
order parameter. Then, it is clear that the order parameter distribution at T « 0.4 shown in 
Fig. undoes not correspond to a critical point of a continuous transition. 

Apparently, the origin of such a peculiar behavior at J2 ~ 0.309 are nietastable states with 
different values of the chiral order parameter which prevent the system investigation considerably. 
The multi-peak structure of the order parameter distribution leads to a sudden jump of the sus- 
ceptibility and states intermediate between mctastable configurations give rise to the specific heat 
anomaly. It should be noted that energy values of the metastable states are close since we observe 
in our simulations that the energy distribution has a Gaussian form even for the largest lattice 
size. 

To estimate the Ising transition temperature at J2 « 0.309 we analyze E'^ distribution defined 
in Eq. ([TO)) and discussed in Sec. 12.3.1 As it is pointed out above, {E'^) ^ whenever a helical 
ordering exists. Fig. [^ shows that the distribution is not symmetric relative to zero at T = 0.52 
and L — 65, while it is definitely symmetric at T > 0.53. Therefore, the critical temperature for 
L = 65 can be roughly estimated as Tj{L) = 0.53(1). The transition temperature can be estimated 
using the relation 

Ti = Tj{L)-^, (26) 

where A is constant and j/ = 1 as it is expected for an Ising transition. An extrapolation to the 
thermodynamic limit using Eq. ([26)) gives Tj = 0.477(12) (see Fig. [23)) . Then, we obtain that 
Tbkt < Ti even in the neighborhood of the Lifshitz point. 



4. Conclusion 

We discuss critical properties of the 2D helimagnet described by the Hamiltonian ([TJ . It belongs at 
J1I J\ > 1/4 to the same class universality as the fully frustrated XY model and the antiferromagnet 
on triangular lattice which have two successive phase transitions upon the temperature decreasing: 
the first one is associated with breaking of the discrete Z2 symmetry and the second one is of the 
BKT type at which the S0{2) symmetry breaks. We confirm that this scenario is realized also in 
the model ([Ij at J2/J1 > 1/4 and obtain the phase diagram shown in Fig. [TJ A narrow region 
exists on this phase diagram between lines of the Ising and the BKT transitions that corresponds 
to the chiral spin liquid. 

In particular, we demonstrate that the number and sequence of transitions do not depend on 
the turn angle ^o of the helix twist at T = 0. Then, this quantity is not a critical parameter 
that has been already found [6l [50] in three-dimensional helimagnets. We find that it is useful 
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in numerical discussion of helimagnets to choose parameters of the Hamihonian so that the hehx 
pitch at T = to be commensurate. It aUows to use definition ^ of the magnetic order parameter 
in which the summation over sublattices is involved. 

We find in accordance with results of Ref. [32] that the specific heat and susceptibilities have 
subsidiary peaks at low T near the Lifshitz point J2/J1 = 1/4 (see Fig. fTT]) . However, in contrast 
to the conclusion of Ref. |32| we demonstrate that these anomalies do not signify a continuous 
phase transition. Apparently, their origin is in metastable states near the Lifshitz point which lead 
also to a peculiar distribution of the chiral order parameter shown in Fig. 1191 

This work was supported by the RF President (grant MD-274.2012.2), the RFBR grant 12-02- 
01234, and the Program "Neutron Research of Solids". 
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